ManyPrimates plots

library(tidyverse)

package ‘tibble’ was built under R version 3.5.2

library(ggthemes)
library(gridExtra)
library(ggpubr) # for ggarrange
library(langcog) # for multi_boot_standard
theme_set(theme_few())
# langcog package is installed by running
# install.packages("devtools")
# devtools::install_github("langcog/langcog")

# make a base plots
# for flipped coordinate system with prop correct on x
p1 <- ggplot() +
  geom_hline(yintercept = 1/3, lty = 2) +
  scale_y_continuous("Proportion Correct", limits = c(-0.05, 1.05)) +
  theme(axis.title.y = element_blank(), axis.ticks.y = element_blank())
p2 <- p1 + theme(axis.text.y = element_blank())

Intro

These are some preliminary visualizations of the data.

phylo <- read_csv("../data/species_data.csv") %>% 
  select(species, species_formatted, clade, clade_formatted, phylo)
mp_data <- read.csv("../data/merged_data/01_manyprimates_pilot_merged_data_v2.csv",
                    stringsAsFactors = F) %>%
  left_join(phylo, by = "species") %>%
  mutate(species = species_formatted,
         species = reorder(species, phylo),
         clade = clade_formatted,
         clade = fct_relevel(clade, "Lemur", "New World monkey", "Old World monkey", "Ape"),
         delay = factor(delay),
         site = fct_recode(site, 
                           "Edinburgh Zoo" = "edinburgh", 
                           "Kumamoto Sanctuary" = "kumamoto sanctuary",
                           "Wolfgang Köhler Primate\nResearch Center" = "Leipzig Zoo",
                           "Sweetwaters Chimpanzee   \nSanctuary Group 1" = "sweetwaters",
                           "Sweetwaters Chimpanzee\nSanctuary Group 2" = "sweetwaters_group2",
                           "Language Research Center" = "LRC"
                           )
)

Overview by species

First an overview of the data, plotted by species and delay. Small, transparent dots represent aggregated data for each individual by delay. Open dots are the group mean for that delay. Error bars are 95% confidence intervals. Text labels are the sample size.

plot_individual <- mp_data %>%
  group_by(phylo, clade, species, delay, subject_site) %>%
  summarise(correct = mean(correct),
            nr_trials = max(trial))
plot_group <- plot_individual %>%
   multi_boot_standard(col = "correct")

|==                                                                                      |  3% ~2 m remaining     
|====                                                                                    |  6% ~48 s remaining    
|==============                                                                          | 17% ~14 s remaining    
|======================                                                                  | 25% ~9 s remaining     
|====================================                                                    | 42% ~4 s remaining     
|===================================================                                     | 58% ~2 s remaining     
|===============================================================                         | 72% ~1 s remaining     
|======================================================================                  | 81% ~1 s remaining     
|================================================================================        | 92% ~0 s remaining     

# get sample sizes (for use in plots)
ns <- mp_data %>%
  group_by(clade, species) %>%
  mutate(n = n_distinct(subject_site)) %>%
  group_by(clade, species, site, n) %>%
  summarise(n_by_site = n_distinct(subject_site))

p2 + aes(x = delay, y = correct, col = delay) +
  geom_jitter(data = plot_individual, aes(size = nr_trials), width = .1, height = .015, alpha = .15) +
  geom_point(data = plot_group, aes(y = mean), shape = 1, size = 3, stroke = 1.5) +
  geom_linerange(data = plot_group, aes(y = NULL, ymin = ci_lower, ymax = ci_upper), lwd = 1.2) +
  # geom_text(data = ns, aes(label = n), y = -.05, x = 2, col = "black", size = 3) +
  facet_grid(species ~ ., switch = "y") +
  theme(strip.text.y = element_text(angle = 180, hjust = 0)) +
  scale_size_area("Number of Trials", max_size = 3, breaks = c(12, 24, 36)) +
  scale_colour_solarized("Delay", breaks = c("short", "medium", "long")) +
  coord_flip(ylim = 0:1)

ggsave("../graphs/02_01_overview.png", width = 4, height = 4.5, scale = 2)
ggsave("../graphs/Fig2.tiff", width = 4, height = 4.5, scale = 2, type = "cairo", compression = "lzw")

Overview by species across delays

plot_individual2 <- plot_individual %>%
  group_by(clade, species, subject_site) %>%
  summarise(correct = mean(correct))
plot_group2 <- plot_individual2 %>%
   multi_boot_standard(col = "correct")

p1 + aes(x = fct_rev(species), y = correct, col = clade) +
  geom_jitter(data = plot_individual2, width = .1, height = .015, alpha = .15, size = 3) +
  geom_point(data = plot_group2, aes(y = mean), shape = 1, size = 3, stroke = 1.5) +
  geom_linerange(data = plot_group2, aes(y = NULL, ymin = ci_lower, ymax = ci_upper), lwd = 1.2) +
  geom_text(data = ns, aes(label = n), y = -.05, col = "black", size = 3) +
  facet_grid(clade ~ ., scales = "free_y", space = "free_y") +
  theme(strip.text.y = element_blank()) +
  scale_colour_solarized() +
  coord_flip()

ggsave("../graphs/02_02_overview_across_delays.png", width = 4, height = 1.8, scale = 2)

Plots by site

Here we select the species for which we have data from multiple sites. This is a very preliminary way of checking whether there is a lot of variation between sites. Plotting conventions are the same as above.

First we check for which species we have data from more than one site:

mp_data %>%
  group_by(species) %>%
  summarise(sites = n_distinct(site)) %>%
  arrange(desc(sites)) %>%
  knitr::kable()

species	sites
Chimpanzee	5
Ring-tailed lemur	2
Brown capuchin monkey	2
Bonobo	2
Gorilla	2
Black-and-white ruffed lemur	1
Black-faced spider monkey	1
Squirrel monkey	1
Rhesus macaque	1
Long-tailed macaque	1
Barbary macaque	1
Orangutan	1

Chimpanzees

chimp_plot_individual <- mp_data %>%
  filter(species == "Chimpanzee") %>%
  group_by(site, delay, subject_site) %>%
  summarise(correct = mean(correct))
chimp_plot_group <- chimp_plot_individual %>%
   multi_boot_standard(col = "correct")

ch <- p2 + aes(x = delay, y = correct, col = delay) +
  geom_jitter(data = chimp_plot_individual, width = .1, height = .015, alpha = .3, size = 3) +
  geom_pointrange(data = chimp_plot_group, aes(y = mean, ymin = ci_lower, ymax = ci_upper), size = .8, shape = 1, stroke = 1.5) +
  geom_text(data = filter(ns, species == "Chimpanzee"), aes(label = n_by_site), y = -.05, x = 2, col = "black", size = 3) +
  facet_grid(site ~ ., switch = "y") +
  theme(strip.text.y = element_text(angle = 180)) +
  scale_colour_solarized("Delay", breaks = c("short", "medium", "long")) +
  ggtitle("Chimpanzees") +
  coord_flip()
ch

Ring-tailed lemurs

rtlemur_plot_individual <- mp_data %>%
  filter(species == "Ring-tailed lemur") %>%
  group_by(site, delay, subject_site) %>%
  summarise(correct = mean(correct))
rtlemur_plot_group <- rtlemur_plot_individual %>%
   multi_boot_standard(col = "correct")

p2 + aes(x = delay, y = correct, col = delay) +
  geom_jitter(data = rtlemur_plot_individual, width = .1, height = .015, alpha = .3, size = 3) +
  geom_pointrange(data = rtlemur_plot_group, aes(y = mean, ymin = ci_lower, ymax = ci_upper), size = .8, shape = 1, stroke = 1.5) +
  geom_text(data = filter(ns, species == "Ring-tailed lemur"), aes(label = n_by_site), y = -.05, x = 2, col = "black", size = 3) +
  facet_grid(site ~ ., switch = "y") +
  theme(strip.text.y = element_text(angle = 180)) +
  scale_colour_solarized(breaks = c("short", "medium", "long")) +
  ggtitle("Ring-tailed lemurs") +
  coord_flip()

ggsave("../graphs/02_04_rtlemur_by_site.png", width = 4, height = 1.2, scale = 2)

Brown Capuchins

cap_plot_individual <- mp_data %>%
  filter(species == "Brown capuchin monkey") %>%
  group_by(site, delay, subject_site) %>%
  summarise(correct = mean(correct))
cap_plot_group <- cap_plot_individual %>%
   multi_boot_standard(col = "correct")

cap <- p2 + aes(x = delay, y = correct, col = delay) +
  geom_jitter(data = cap_plot_individual, width = .1, height = .015, alpha = .3, size = 3) +
  geom_pointrange(data = cap_plot_group, aes(y = mean, ymin = ci_lower, ymax = ci_upper), size = .8, shape = 1, stroke = 1.5) +
  geom_text(data = filter(ns, species == "Brown capuchin monkey"), aes(label = n_by_site), y = -.05, x = 2, col = "black", size = 3) +
  facet_grid(site ~ ., switch = "y") +
  theme(strip.text.y = element_text(angle = 180)) +
  scale_colour_solarized("Delay", breaks = c("short", "medium", "long")) +
  ggtitle("Capuchin Monkeys") +
  coord_flip()

lay <- matrix(c(1, 1, 2))

grid.arrange(ch, cap, layout_matrix = lay)

ggsave("../graphs/02_05_chimp_capuchin_by_site.png", arrangeGrob(ch, cap, layout_matrix = lay), width = 4, height = 4, scale = 2)
ggsave("../graphs/Fig5.tiff", arrangeGrob(ch, cap, layout_matrix = lay), width = 4, height = 4, scale = 2, type = "cairo", compression = "lzw")

Bonobos

bon_plot_individual <- mp_data %>%
  filter(species == "Bonobo") %>%
  group_by(site, delay, subject_site) %>%
  summarise(correct = mean(correct))
bon_plot_group <- bon_plot_individual %>%
   multi_boot_standard(col = "correct")

p2 + aes(x = delay, y = correct, col = delay) +
  geom_jitter(data = bon_plot_individual, width = .1, height = .015, alpha = .3, size = 3) +
  geom_pointrange(data = bon_plot_group, aes(y = mean, ymin = ci_lower, ymax = ci_upper), size = .8, shape = 1, stroke = 1.5) +
  geom_text(data = filter(ns, species == "Bonobo"), aes(label = n_by_site), y = -.05, x = 2, col = "black", size = 3) +
  facet_grid(site ~ ., switch = "y") +
  theme(strip.text.y = element_text(angle = 180)) +
  scale_colour_solarized(breaks = c("short", "medium", "long")) +
  ggtitle("Bonobos") +
  coord_flip()

ggsave("../graphs/02_06_bonobo_by_site.png", width = 4, height = 1.2, scale = 2)

Gorilla

gor_plot_individual <- mp_data %>%
  filter(species == "Gorilla") %>%
  group_by(site, delay, subject_site) %>%
  summarise(correct = mean(correct))
gor_plot_group <- gor_plot_individual %>%
   multi_boot_standard(col = "correct")

p2 + aes(x = delay, y = correct, col = delay) +
  geom_jitter(data = gor_plot_individual, width = .1, height = .015, alpha = .3, size = 3) +
  geom_pointrange(data = gor_plot_group, aes(y = mean, ymin = ci_lower, ymax = ci_upper), size = .8, shape = 1, stroke = 1.5) +
  geom_text(data = filter(ns, species == "Gorilla"), aes(label = n_by_site), y = -.05, x = 2, col = "black", size = 3) +
  facet_grid(site ~ ., switch = "y") +
  theme(strip.text.y = element_text(angle = 180)) +
  scale_colour_solarized(breaks = c("short", "medium", "long")) +
  ggtitle("Gorillas") +
  coord_flip()

ggsave("../graphs/02_07_gorilla_by_site.png", width = 4, height = 1.2, scale = 2)

Task-experience

Here we split each species by task experience. Check if we have species with sufficient memebers having different levels of task experience.

mp_data %>%
  group_by(species) %>%
  mutate(lvls_task_exp = n_distinct(task_experience)) %>%
  filter(lvls_task_exp > 1) %>%
  group_by(species, task_experience) %>%
  summarise(n = n_distinct(subject_site)) %>%
  knitr::kable()

species	task_experience	n
Chimpanzee	no	19
Chimpanzee	yes	32
Bonobo	no	6
Bonobo	yes	5

# get sample sizes (for use in plots)
ns_task_exp <- mp_data %>%
  group_by(clade, species) %>%
  mutate(n = n_distinct(task_experience)) %>%
  group_by(clade, species, task_experience, n) %>%
  summarise(n_by_site = n_distinct(subject_site))

So far, this only applies to chimps and bonobos. However, task experience co-varies with site.

chimp_task_plot_individual <- mp_data %>%
  filter(species == "Chimpanzee") %>%
  group_by(site, task_experience, delay, subject_site) %>%
  summarise(correct = mean(correct))
chimp_task_plot_group <- mp_data %>%
  filter(species == "Chimpanzee") %>%
  group_by(task_experience, delay, subject_site) %>%
  summarise(correct = mean(correct)) %>%
  multi_boot_standard(col = "correct")

p_taskexp_chimp <- p2 + aes(x = delay, y = correct) +
  geom_jitter(data = chimp_task_plot_individual, aes(col = site), width = .3, height = .015, alpha = .5, size = 3) +
  geom_pointrange(data = chimp_task_plot_group, aes(y = mean, ymin = ci_lower, ymax = ci_upper, shape = delay), size = .8, stroke = 1.5) +
  geom_text(data = filter(ns_task_exp, species == "Chimpanzee"), aes(label = n_by_site), y = -.05, x = 2, col = "black", size = 3) +
  facet_grid(task_experience ~ ., switch = "y") +
  theme(legend.box = "horizontal", strip.text.y = element_text(angle = 180)) +
  scale_shape_manual(values = c(1, 2, 5), breaks = c("short", "medium", "long")) +
  scale_colour_solarized() +
  ggtitle("Chimpanzees") +
  coord_flip()

bonobo_task_plot_individual <- mp_data %>%
  filter(species == "Bonobo") %>%
  group_by(site, task_experience, delay, subject_site) %>%
  summarise(correct = mean(correct))
bonobo_task_plot_group <- mp_data %>%
  filter(species == "Bonobo") %>%
  group_by(task_experience, delay, subject_site) %>%
  summarise(correct = mean(correct)) %>%
  multi_boot_standard(col = "correct")

p_taskexp_bon <- p2 + aes(x = delay, y = correct) +
  geom_jitter(data = bonobo_task_plot_individual, aes(col = site), width = .3, height = .015, alpha = .5, size = 3) +
  geom_pointrange(data = bonobo_task_plot_group, aes(y = mean, ymin = ci_lower, ymax = ci_upper, shape = delay), size = .8, stroke = 1.5) +
  geom_text(data = filter(ns_task_exp, species == "Bbonobo"), aes(label = n_by_site), y = -.05, x = 2, col = "black", size = 3) +
  facet_grid(task_experience ~ ., switch = "y") +
  theme(legend.box = "horizontal", strip.text.y = element_text(angle = 180)) +
  scale_shape_manual(values = c(1, 2, 5), breaks = c("short", "medium", "long")) +
  scale_colour_solarized() +
  ggtitle("Bonobos") +
  coord_flip()

grid.arrange(p_taskexp_chimp, p_taskexp_bon, ncol = 1)

grob <- arrangeGrob(p_taskexp_chimp, p_taskexp_bon, ncol = 1)
ggsave("../graphs/02_08_task_experience.png", grob, width = 4, height = 3, scale = 2)

Age

Here we plot age against correct choice separate for each delay and species. Regression line is smoothed delayal mean.

plot_age <- mp_data %>%
  mutate(delay = fct_rev(delay)) %>%
  group_by(subject_site, norm_age, clade, species, delay) %>%
  summarise(correct = mean(correct))

ggplot(plot_age, aes(x = norm_age, y = correct)) +
  geom_jitter(aes(col = clade), width = .05, height = .05, alpha = .5, size = 2.5) +
  geom_smooth(method = "glm", method.args = list(family = "binomial"), col = "black") +
  # geom_vline(xintercept = 0, lty = 2) +
  geom_hline(yintercept = 1/3, lty = 2) +
  facet_grid(~ delay) +
  labs(x = "Normed Age (relative to species longevity)", y = "Proportion Correct") +
  scale_color_solarized("Clade") +
  ylim(c(-.05, 1.05))

ggsave("../graphs/02_09_age_by_delay.png", width = 4, height = 1.3, scale = 2)
ggsave("../graphs/Fig4.tiff", width = 4, height = 1.3, scale = 2, type = "cairo", compression = "lzw")

ggplot(plot_age, aes(x = norm_age, y = correct)) +
  geom_jitter(aes(fill = species), width = .05, height = .05, alpha = .5, size = 2.5, shape = 21, stroke = 0) +
  geom_smooth(aes(col = delay), method = "glm", method.args = list(family = "binomial"), show.legend = F) +
  geom_vline(xintercept = 0, lty = 2) +
  geom_hline(yintercept = 1/3, lty = 2) +
  facet_grid(delay ~ clade, scales = "free_x") +
  labs(x = "Normed Age (relative to species longevity)", y = "Proportion Correct") +
  theme(legend.position = "bottom") +
  scale_color_manual(values = rev(solarized_palette(3))) +
  ylim(c(-.05, 1.05))

ggsave("../graphs/02_10_age_by_delay_species.png", width = 4, height = 3.4, scale = 2)

Cup distance and board size

(sig. effect in model)

plot_cup <- mp_data %>%
   group_by(subject_site, cup_distance, clade, species, delay) %>%
   summarise(correct = mean(correct))

cp <- ggplot(plot_cup, aes(x = cup_distance, y = correct)) +
  geom_jitter(aes(col = clade), width = .5, height = .05, alpha = .5, size = 2) +
  geom_smooth(method = "glm", method.args = list(family = "binomial"), col = "black") +
  labs(x = "Cup Distance in cm", y = "Proportion Correct") +
  geom_hline(yintercept = 1/3, lty = 2) +
  scale_color_solarized("Clade") +
  ylim(c(-.05, 1.05)) +
  ggtitle("A. Cup Distance")

plot_board <- mp_data %>%
   group_by(subject_site, board_size, clade, species, delay) %>%
   summarise(correct = mean(correct))

bp <- ggplot(plot_board, aes(x = board_size, y = correct)) +
  geom_jitter(aes(col = clade), width = .5, height = .05, alpha = .5, size = 2) +
  geom_smooth(method = "glm", method.args = list(family = "binomial"), col = "black") +
  labs(x = "Board Size in cm", y = "Proportion Correct") +
  geom_hline(yintercept = 1/3, lty = 2) +
  scale_color_solarized("Clade") +
  ylim(c(-.05, 1.05)) +
  ggtitle("B. Board Size")

pcb <- ggarrange(cp, bp, common.legend = T, legend = "right")

non-integer #successes in a binomial glm!non-integer #successes in a binomial glm!non-integer #successes in a binomial glm!

pcb

ggsave("../graphs/02_11_cup_distance_board_size.png", pcb, width = 6, height = 2.3, scale = 2)
ggsave("../graphs/FigS1.tiff", pcb, width = 6, height = 2.3, scale = 2, type = "cairo", compression = "lzw")